%%

% try different shifts
% include the REST class to soak variance 
   % - in the movement / non movement check

%% segement all the files
FileList = {
   
  %  'mohEEG111.mat'
   % 'Marion12_2504.mat'
   % 'Tim12_2804.mat'
   %   'Romain12_2904.mat'
      'david12-2904.mat'
    };

%FileList = { 'session2.acq.mat' };
%Fs = 2000; default
Fs = 2000;

window = 520; % 3000 = 1.5 seconds
window_buffer= 5500; % not used 
shift = 0;
X = []; % contains the trials from all classes
Y = []; % contains the classification from all classes

min_window = 800; % samples
overlap = window-window/8;
d_min = ceil ( ((min_window - window) / overlap) + 1 );

for curFileIx = 1:length(FileList)
     
        load(['D:\temp\subjects\' FileList{curFileIx}]);
        disp(['Subject file:' FileList{curFileIx}]);

        datar = data(:,1:size(data,2)-2); % remove 2 trigger channel at the end

        %trial_channel = data(:,8);
        trial_channel = data(:, size(data,2)-1);

        %class_separation_channel = data(:,9);
        class_separation_channel = data(:, size(data,2));

        classN = length(count_sequence(class_separation_channel,5));  % "5" tags each class


        classRest = classN; % the rest class is always the last

        trialsN = length(count_sequence(trial_channel,5));  % "5" tags each trial
            
        disp(['Classes/Movements detected: ' int2str(classN)]);
        disp(['Total trails detected: ' int2str(trialsN)]);
        disp(['Frequency: ' int2str(Fs)]);
        
        for m = 1:classN 
            
            [best_index,best_trigger,class_start_pos,class_end_pos] = generate_trials_postitions_using_discrete_channel(m, trial_channel,class_separation_channel);
            positions{m} = [class_start_pos,class_end_pos];
            
            if (m<classN) %we want to take the signal only for the rest
                continue;
            end;
            
            % cut signal
            %index = find(diff(trigger)==1); % because we do it earlier now!!!!
            trials_per_class = zeros(size(datar,2),window,length(best_index));
            disp(['Trials: ' int2str(length(best_index)) ' ' FileList{curFileIx}])

            %get the data for these triggers
            for i=1:length(best_index)
                %calculate a range around the trigger
                buffer = datar(best_index(i)+shift:best_index(i)+shift+window-1,:)';
                trials_per_class(:,:,i) = buffer;
            end

            X = cat(3,X,trials_per_class); % put all trials together
            Y = cat(1,Y,ones(length(best_index),1) * m);

        end;
   
       disp('--------------------------------');% next file
end

%% detect movements
figure;

for iclass = 1:11
    
    subplot(3,4,iclass);
    
    firstClass = datar(positions{iclass}(1):positions{iclass}(2),:);

    %% we train a potato to detect movements

    COV_rest=eeg2cov(X,window,overlap);

    [C, th] = potato_estimation_iterativ(COV_rest);

    [COVtest CovTestIndex] = eeg2cov_index(firstClass',window,overlap);

    [~,d] = potato_detection(COVtest,C,th);

%    t = size(rms(firstClass(:,1), window, overlap, 0),2);
%    data_average = zeros(t,size(firstClass,2));

%     for i = 1:size(firstClass,2)
%         averaged = rms(firstClass(:,i), window, overlap, 0);
%         data_average(:,i) = averaged;
%     end;

    %plot(data_average);

%    disp('Done.')

    %% find good threshold 
    
    xx = (1:length(d))';
    f=fit(xx,d,'poly2');

    f = f(xx);
    f = f * (9/10); % correction of threshold
    above_d = d > f;

    [seq pos]= count_sequence(above_d,1); % position in ds

    m = mean(seq);
    st = std(seq);
    minn1 = m - 1 * 2 * st;
    minn2 = m + 2 * st;

    temp = 0;
    movements = 0;

    removed = 0;

    limit = 30;
    
    %% remove bad movements 
    
    good_ones = seq > limit; %seq > minn1 & seq < minn2;
    seq_good = seq(good_ones);
    pos_good = pos(good_ones);

    bad_ones = (seq <= limit);
    seq_bad = seq(bad_ones);
    pos_bad = pos(bad_ones);

    d_fixed = d;
    
    %actual removal
    for i = 1:length(seq_bad)

        d_fixed(pos_bad(i):pos_bad(i) + seq_bad(i)) = 1;
    end;

    result_binary = zeros(positions{iclass}(2)-positions{iclass}(1),1);

    %% extract signal for each movement
    
    for i = 1:length(seq_good)

        p = pos_good(i) + seq_good(i) / 3; % find starting position in d space

        k = CovTestIndex(pos_good(i)); % find the exact position in sample space 

        window_buffer = ceil( seq_good(i) * 53 * (4/10) ); % reduce by 10 percent;

        buffer = firstClass(k:k+window_buffer,:);
        result_binary = result_binary + [zeros(k,1);ones(window_buffer,1);zeros(length(result_binary) - (k+window_buffer),1)];
        movements = movements + 1;
    end;
    
    movements
    
    %% Draw the distances 
    %  figure;
    % plot(result_with_min);

%     plot(xx,d_fixed(1:length(xx)));
%     hold on;
%     plot(xx,d_fixed(1:length(xx)) .* (d_fixed(1:length(xx))>f(xx)),'r');
    
%    figure;
    plot(trial_channel(positions{iclass}(1):positions{iclass}(2),:),'g');
    hold on;
    plot(result_binary*2);
    ylim([-1,7]);
    
end; % end of class loop



